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A STUDY OF THE ACCURACY OF ESTIMATING THE ORBITAL 


ELEMENTS OF A LUNAR SATELLITE BY USING RANGE 

AND RANGE -RATE MEASUREMENTS 

By Harold R. Compton 
Langley Research Center 


SUMMARY 


A parametric study has been made to determine the effects of tracking mode 
and orbital parameters on the accuracy of determining the state of a lunar sat- 
ellite. A comparison of the relative advantage of using either range or range- 
rate measurements was made^ and the results indicate that the advantage of one 
data type over the other is very dependent on the semimajor axis and somewhat 
dependent on the nodal position and eccentricity but is not dependent on the 
inclination. For example^ it was found that for lunar orbits with medium eccen 
tricity and with semimajor axes of approximately 2500 kilometers, the orbital 
elements could be determined with equal accuracy by using either range or range 
rate data when the ratio of the standard deviation of range measurements to the 
standard deviation of range-rate measiirements was approximately I 5 OO seconds. 
For lunar orbits with semimajor axes smaller than approximately 2500 kilometers 
range-rate measurements gave a more accurate determination of the elements, 
whereas for orbits with semimajor axes greater than approximately 2500 kilo- 
meters, range measurements proved to be the better data type of the two. It 
was also found that range and range -rate measiirements are similar data types in 
that they produce similar correlation matrices for the state variables and 
simultaneous use of both data types does not significantly reduce correlations 
between the elements. 

Over a range of inclination angles from 2^ to 40^ with respect to the 
earth-moon plane, the results showed that the accuracy of determining the ori- 
entation angles increased as the inclination increased whereas the accuracy of 
determining the in-plane variables remained nearly constant. The effects of 
variations in the nodal position on the accuracy of determining the orbital 
elements were found to be periodic. The eccentricity, argument of periapsis, 
and the longitude of the ascending node were best determined when the orbit was 
viewed on edge whereas the semimajor axis, inclination, and the time of periap- 
sis passage were best determined when the orbit was viewed broadside. It was 
also found that the accuracy of estimating the elements increased with an 
increase in eccentricity. 



IliTTRODaCTION 


Current plans for lunar research missions include the establishment of 
satellites in orbit about the moon. For unmanned missions the elements of the 
orbit in which the vehicle is moving must be known within a reasonable degree 
of accuracy in order to determine the location of the satellite when data are 
taken by the satellite. For manned missions the position must be known with a 
high degree of accuracy. The knowledge of the orbital elements may also be 
used to determine certain selenodetic constants and in particular the coeffi- 
cients of the harmonics of the lunar gravitational field. It is therefore of 
interest to investigate the accuracy to which the orbital elements of a l\onar 
satellite can be determined by earth-based tracking. 

The basic earth -based data types are range ^ range-rate^ and angular meas- 
urements. Since the accuracy of making angular measurements of a vehicle 
moving in orbit about the moon is low^ this particular data type was not con- 
sidered. Hence the results presented in this paper are based solely on the use 
of range and range-rate measurements. It is of interest not only to estimate 
the accuracy to which the elements can be determined but also to ascertain the 
relative advantage of one data type over the other. Therefore a parametric 
study in which both data types were used was initiated. 

In order to make the parametric study^ the statistical equations which 
were used to estimate the accuracy to which the orbital elements could be deter- 
mined were programed in double precision on an IBM 709^ electronic data proc- 
essing system. The basic theory and equations used in the program are given in 
appendixes A and B. 


SYMBOLS 


Unless otherwise specified^ the unit of length is the lunar radius^ which 
is 1738 kilometers^ and the unit of time is the period of a lunar surface 
satellite divided by which is 1035 seconds. The coordinate system and 
angular parameters are illustrated in figure 1 . 

A,B matrices containing partial derivatives of a given data type with 

respect to orbital elements 

a semimajor axis of lunar satellite orbit 

E eccentric anomaly; operator used in appendix A 

e eccentricity of lunar satellite orbit 

fj_ functional relation between observable quantity and parameters to 

be estimated (see eq. (A 2 )) 

i inclination of orbital plane of lunar satellite to earth-moon plane 
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direction cosines (see appendix B) 

mean anomaly 

munber of observations 

mean angular rate of lunar satellite 

position of lunar satellite 

weighted least squares function defined in equation (AJ) 
distance from center of earth to center of moon 
distance from center of moon to lunar satellite 
vector from center of moon to lunar satellite 
time 

time of ith measurement 
time of periapsis passage 
true anomaly 
weighting matrix 

coordinate axes with origin at center of moon (The X-axis is posi- 
tive in the direction from the center of the moon away from the 
center of the earthy the Y-axis is positive in the direction of 
rotation of the moon^ and the Z-axis is positive in such a direc- 
tion that it forms a right-handed axis system.) 

position components of lunar satellite 

ith measurement of general quantity y, where i = 1^ 2, . . . N 
parameter to be estimated (subscript denotes particular parameter) 
error 

error in ith measurement 
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Aq, covariaxice matrix of estimated parameters 

gravitational constant of moon 

I mean motion of moon about earth 

p range or distance from center of earth to position of lunar sat- 

ellite j in appendix A the symbol p is used with double sub- 
scripts to denote the correlation between the variables indi- 
cated by the subscripts 

p^ ith measurement of p 

■p vector from center of earth to position of lunar satellite 

p range rate or radial velocity of lunar satellite with respect to 

center of earth 

p^ ith measurement of p 

a standard deviation or one-sigma uncertainty (When this symbol 

appears with a subscript, it is taken to mean the one-sigma 
uncertainty in the estimation of the variable indicated by the 
subscript. ) 

CD argument of perlapsis, angle measured in luna?* satellite plane 

from ascending node to perlapsis 

longitude of ascending node of lunar satellite orbital plane 
measured in earth -moon plane in direction of rotation of moon 
from positive X-axis 

n' longitude of ascending node meas\ired in YZ-plane (see sketch l) 


AlIALYSIS 


In order to simplify the problems associated with the analysis in this 
study, certain assumptions were made. The moon was assumed to be a point mass 
rotating about the earth in a circular orbit. A single observation station 
making uncorrelated, unbiased range and range-rate measurements of a lunar sat- 
ellite moving in a two-body orbit and not occulted by the moon was assumed to 
be located at the center of the earth. All the results presented were obtained 
by assuming a constant one-sigma error in the range measurements and range-rate 
measurements of I5 meters and 0.01 meter per second, respectively. These val- 
ues are conservative estimates of the tracking data accuracy applicable to the 
NASA deep space net (DSN) tracking system (see ref. l). Data were ass\imed to 
have been processed for 1, 2, 4, sometimes 10 consecutive orbits. 

Usually 26 range and 26 range-rate measurements equally spaced in time were 
simulated during each orbit. 
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The fundamental approach made in this error analysis was to simulate range 
and range-rate measurements over a given period of time and from this simulation 
to calculate a covariance matrix from which the variances of the elements could 
he obtained. The covariance matrix was obtained from a weighted least squares 
simulation. (See appendix A for details of the simulation.) The form of this 
matrix is 


Ax = a2(ATA)"^ (l) 

where the (i^j) element of A is equal to the partial derivative of the ith 
observation with respect to the jth element to be estimated and where a is 
the standard deviation of the measurements. Explicit expressions for the par- 
tials of range and range rate with respect to the Keplerian elements are given 
in appendix B. It can be seen from equation (l) that no actual values of the 
measurements are needed^ and in particular only the standard deviation of the 
observations is used. Therefore, with a fixed tracking schedule and only one 
type of data, the one-sigma uncertainty in the estimation of the elements is 
proportional to the standard deviation of the data type. Thus, when the accu- 
racy of estimating the elements is compared for two data types, an important 
parameter is the ratio of the standard deviations of the data types, for 
exanple, The covariance matrix for the simultaneous use of two data 

types is shown in appendix A. 

In order to make a parametric study of the effects of a given element on 
the accuracy of determining the elements, five elements of a chosen nominal 
orbit were held constant and the sixth was varied over a given range. The one 
exception to this procedure was in the eccentricity variation for which, instead 
of the nominal value for the semimajor axis, a value of 5000 kilometers was 
used. This exception was made in order to insure that over the given range of 
eccentricities the distance from the center of the moon to the lunar satellite 
was never less than the radius of the moon. The elements of the nominal orbit 
used in this investigation were chosen to provide a low periapsis (approximately 
50 km) in an orbit with medium eccentricity and inclination. These elements are 
as follows: 


a = 2255 kilometers 
i - 30 ^ 
a = 300 

CO = l80^ 

e = 0.2 

to = 0 second 


In this study the angle co always appears in the partial derivatives as 
an angle added to the true anomaly v in the argument of either a sine or 
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cosine function. Inasmuch as the true anomaly rotates through 560^ each orbit, 
the argument of this sine or cosine function rotates through one period regard- 
less of the value of o). Hence, the effects of a variation in o) upon the 
accuracy of estimating the elements are negligible as long as integral orbits 
of tracking are used and occultations are not considered. Similarly, no results 
are presented for a variation in to^ because this parameter only defines vhere 
the vehicle is located in the orbit and, 'as long as measurements are made over 
complete orbits, the effects of changing tQ are negligible. 

RESULTS AM) DISCUSSION 


Effects of Tracking Schedule on Accuracy of Estimating the Elements 

The partial derivatives contained in the A matrix are fundamental to the 
entire orbit -determination process. A large derivative is said to have a large 
information content and, similarly, a small derivative is said to have a small 
information content. It can be shown that the accioracy of estimating any param- 
eter increases as the information about that parameter increases or hence as the 
derivative of the observable quantity (range or range rate) with respect to the 
parameter increases. Figures 2 and 5 are presented to show how the partials of 
range and range rate with respect to the orbital elements change with time. 

These derivatives are plotted as functions of time over five orbital periods. 

The orbital period is approximately 2.6 hours. It can be seen that the partials 
have a periodic nature with a period equal to that of one orbit. 

In figures 2 and 5 the amplitudes of the curves representing the partial 
derivatives of range and range rate with respect to the semimajor axis continue 
to increase with time. This increase is due to the mixed secizlar terms such as 
n(t - to) sin E which are contained in the analytical expressions for these 
derivatives. Because the information content of the data increases rapidly 
with time due to the mixed secular terms, the semimajor axis shoiild be deter- 
mined more accurately over long time arcs. The amplitudes of the partial deriv- 
atives of range and range rate with respect to the inclination appear to be 
decreasing with time, but the reason for this decrease is the fact that the 
derivatives vary as sin Q and, in the particular case shown in figures 2 
and 5^ ^ is 30^ at time zero and is decreasing at a rate of 0.5^® per hour 
due to the rotation of the X-axis. Since the amplitude of the partial deriva- 
tive with respect to inclination is a maximum when 0 , = 90^ or = 27OO, it 
is expected that the inclination would be best deteimined when the orbit is 
viewed broadside, that is, when Q equals 90^ or 270^. Likewise, a very weak 
determination of the inclination is expected when 2 = 0^ or 0 , - l80^ 
because, regardless of the inclination, identical time histories of range and 
range-rate measurements would be obtained - that is, the observations are inde- 
pendent of inclination. 

When the orbit is viewed nearly on edge, that is, when Q is very near 0^ 
or 180^, the partial derivatives of range and range rate with respect to Cl 
and 0) are related through the expressions 
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Therefore if the inclination is not large ^ the derivatives are approximately 
equal, and figures 2 and 3 show that even with an inclination of ^0^ and with 
Q = 30^ they are not vastly different. If the relations in equations (2) 
and ( 3 ) were always exactly true, the normal matrix A^a/ct^ would he singular 
and noninvertihle , because one column of the A matrix would he proportional 
to another column and therefore one row of the normal matrix would he propor- 
tional to another row. Thus, in a real-orhit -determination process where the 
measurements are made with the orhit being viewed nearly on edge, the noimial 
matrix is expected to he poorly conditioned for inversion by use of finite- 
decimal arithmetic. 

The special case of near-zero inclination also has problems associated 
with it. By referring to the equations in appendix B, it can he seen that when 
i is nearly zero the partial derivatives of range and range rate with respect 
to n and O) are related through the following expressions regardless of the 
nodal position: 




( 4 ) 





( 5 ) 




Again, the normal matrix is expected to he poorly conditioned for inversion. 
Thus, trouble might he expected in trying to invert the normal matrix asso- 
ciated with an orhit having near-zero inclination. 

No physical significance should he attached to the fact that the partial 
derivative of range with respect to the eccentricity remains positive as shown 
in figure 2. This fact is due to the position of the line of nodes during the 
observation period. The node angle would he different for different observa- 
tion periods and hence the derivative might he negative. 

The accuracy of estimating the orbital elements varies significantly with 
the number of orbits tracked, the number of observations made during each orbit, 
and the position of the satellite in the orbit at the time of the observations. 
In order to show how the accuracies vary with the number of orbits tracked, 
figure k is presented. This figure illustrates how well the elements of the 
nominal orbit could have been determined if the vehicle had been tracked by 
making range and range-rate measiirements every 6 minutes over a period of 1 to 
10 orbits. As would be expected, the accuracy of estimating the elements 
increases with the number of orbits tracked. It is also seen in figure k that 



range measurements and range-rate measurements having one-sigma errors of 
15 meters and 0.01 meter per second^ respectively, (c^p/o'p = 1500 seconds) can 

be used separately to determine the elements to approximately the same accuracy. 
It is shovn subsequently that this result is due to the particular values of 
the orbital elements utilized and that in particular the relative advantage of 
one data type over the other is dependent upon the semimajor axis, the nodal 
position, and eccentricity. Thus, this result cannot be generalized to include 
all lunar orbits. The lowest cirrve in each set of curves presented in figure k 
shows the one- sigma values resulting from the simultaneous use of range and 
range-rate measurements. These values are, of course, smaller than the values 
obtained when either data type is used separately. 

Interpolation in figure ^4- between integral orbits to obtain one-sigma val- 
ues for fractions of orbits is only approximate, and the reason for this can be 
understood by referring to the discussion of figures 2 and 5 where it was shown 
that, in general, the time variation of the partial derivatives is periodic, 
with a period being equal to that of one orbit. Since the elements of the nor- 
mal matrix are sums of products of these partial derivatives and since the 
covariance matrix is the inverse of the normal matrix, it is apparent that the 
values of the elements of the covariance matrix are largely dependent upon the 
times at which the observations were made. Therefore if measurements are simu- 
lated over a fraction of an orbit, the values of the matrix elements will vary 
according to the portion of the orbit investigated, and for this reason the 
curves in figure ^ between integral orbits may not be as smooth as indicated. 
Hence, interpolation in figure 4 is only approximate. 

Before proceeding to other results, it should be noted that in a real orbit 
determination the covariance matrix associated with tracking the satellite for 
a single orbital period may be difficult to obtain due to mxmerical operations 
in the con5)uter. In particular, the normal matrix may be nearly singular, and 
experience has shown that single -precision 
arithmetic is not adequate for inversion 2. 

of this matrix. During a single orbital 
period, the moon rotates through a very 
small angle. Thus, the problem under con- 
sideration approaches the stationary -moon 
problem for which the normal matrix A^a/ot^ 
becomes singular. The covariance matrix is 
the inverse of the normal matrix, as was 
shown in equation (l), and hence it cannot 
be obtained if the normal matrix is 
singular . 

The following argument shows why the 
normal matrix is singular when the moon is 
stationary. It has been shown in refer- 
ence 2 that, in the case of the stationary 
moon, range -rate measurements (and it can 
be similarly shown for range measurements) X 
are independent of the angle which is 

defined in sketch 1. Sketch 1 
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Therefore, the following relations hold: 


= 0 

(6) 

5a' 


=0 

(7) 




By writing i, Q, and (d as functions of ^2’ and noting that p and p 
are functions of i, ^2, and o), the following equations can be derived from 
eq\iations (6) and ( 7 ) by direct substitution: 


cos n 



COS i sin H Sp ^ 

sin n 8p __ Q 

(8) 


sin i 

sin i ^ 


cos i sin ^2 8p 

sin n Sp ^ 

(9) 


sin i 

sin i ScD 


The coefficients of the partial derivatives in equations ( 8 ) and ( 9 ) are the 
partial derivatives of i, and co with respect to $2^ . Equation ( 8 ) 

indicates that 8 p/ 8 i^ 8 p/ 8 f 2 , and Bp/3oD are linearly related and hence 
the A matrix is at most of rank 5- Thus, the normal matrix A^a/u^ has 
a rank of at most 5 and is therefore singular and noninvertible. From 
eq\oation ( 9 ) it is clear that the same conclusions hold for range-rate data. 
Usually if a normal matrix is poorly conditioned for inversion, high corre- 
lations between the parameters which cause this poor conditioning can be 
expected. A high correlation between i, and a> is therefore expected, 

and, as shown subsequently, these three parameters are highly correlated. 

The problem of detemnining the elements based on range and range-rate data 
for a single orbital period is not exactly the stationary-moon problem because 
the moon has rotated through an angle of 1.5^. However, it is questionable 
whether 1.5^ of rotation is sufficient to reduce the linear relations between 
i, n, and cd enough to allow the normal matrix to become invertible in finite- 
decimal arithmetic. Experience has shown that single-precision, 8-decimal 
arithmetic is not adequate to Invert the normal matrix associated with one 
orbit of tracking. This difficulty was circumvented by programing the problem 
in double precision, that is, l6-decimal arithmetic. 

It is of interest to know whether the one-sigma estimation error presented 
in this report obeys the \fn law which states that the one-sigma uncertainty 
in estimating a parameter is inversely proportional to the ^ where N is 
the number of measurements made. By referring to equation (l) , it can be seen 
that if the partial derivatives which are used in the elements of A*^A were 
constant, then the one-sigma estimation error would obey the law exactly. 

It has been shown that these derivatives are not constant but periodic and 
therefore if a sufficient number of measurements were used over the tracking 
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interval, one vould expect the elements of A^A to be approximately propor- 
tional to N. Hence one would expect the one-sigma estimation error to approx- 
imately obey the >/N law. Figure 5 is presented to show that the data in this 
report do approximate this law and hence can be generalized to include an arbi- 
traiy number of observations. In order to obtain the one-sigma values shown in 
figure 5^ was assumed that N observations of the satellite were made over 
five consecutive orbits. This process was repeated several times for the same 
five orbits but with N changed each time. These one-sigma values were 
plotted as f\anctions of N on log-log paper in figure 5^ and it can be seen 
that for all six elements the curves have a slope of approximately -I/2, as 
was expected. 

In the discussion of figure h it was pointed out that the simultaneous 
use of the two data types produced estimates of the elements which were more 
accurate than those obtained from the use of either data type alone. However, 
by referring to figure 5 it can be seen that an accuracy equivalent to that 
obtained from the simultaneous use of the two data types can be achieved by 
using more observations of the same data type. 

In a real orbit determination, the simultaneous use of several data types 
would be expected to help eliminate high correlations between the parameters 
and thereby cause the normal matrix to be better conditioned for inversion. 

As stated previously, if a normal matrix is poorly conditioned for inversion, 
high correlations between the parameters causing the poor conditioning can 
usually be expected. Therefore, it is of interest to know whether the simul- 
taneous use of range and range-rate data serves to eliminate high correlations 
between the parameters. 

The correlation matrices obtained after one orbit by using range, range- 
rate, and range plus range-rate data are presented in figure 6 (see appendix A 
for the definition of the correlation matrix). In the previous discussion of 
the stationary -moon problem, it was pointed out that i, Qy and O) are 
expected to be highly correlated after one orbit of tracking, and figirre 6 
shows that a high correlation does exist. Note that the correlations between 
these three parameters which were obtained by the simultaneous use of range 
and range-rate data are not significantly lower than those obtained when either 
data type is used alone. The correlation matrices after five orbits are shown 
in figure 7. It can be seen that the parameters i, H, and o> are still 
highly correlated in all three correlation matrices. The correlation between 
0 , and o) is slightly higher than that between i and Q and i and co 
due to the nearly linear relation between Q and co when i is small. A 
comparison of the three correlation matrices was made after each orbit up to 
10 orbits. Except for the one-orbit case, it was found that the three matrices 
were similar - that is, elements which were highly correlated on one matrix 
were highly correlated on the other. Hence it was concluded that range and 
range-rate measurements are similar data types and that their simultaneous use 
does not produce any appreciable reduction in the correlations. 
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Effects of Variations in the Elements on the 


Accuracy of Determining the State 

The effects upon the accuracy of determining the elements due to a varia- 
tion in the semimajor axis are illustrated in figure 8 , This figure is a plot 
of the one-sigma value of the error in the estimation of the elements after 
five consecutive orhits of tracking as a function of the semimajor axis. The 
semimajor axis was varied between 2255 5000 kilometers^ and i^ Qy e, 

and to were held constant at the nominal values previously given. For each 
value of the semimajor axis considered in the present study ( 2235 ^ 2500 ^ 3500 ^ 
4000^ 4500 ^ and 5 OOO km) ^ it was assumed that 26 range and 26 range-rate meas- 
urements equally spaced in time were made every orbit for five consecutive 
orbits. It can be seen that^ except for the elements a and to^ the accuracy 
of determining the elements increased as the semimajor axis increased when 
range measurements alone were used^ whereas the accuracy decreased for the 
elements a^ e^ and to and increased for the elements i^ and 03 when 

range-rate measurements alone were used. In general^ for lunar orbits with 
semimajor axes of approximately 2500 kilometers^ the orbital elements could be 
determined with equal accuracy by using either range or range -rate data when 
cfp/a^ ^ 1500 seconds. For lunar orbits with semimajor axes smaller than 

approximately 25 OO kilometers^ range-rate measiirements gave a more accurate 
determination of the elements, whereas for orbits with semimajor axes greater 
than approximately 2500 kilometers, range measurements proved to be the better 
data type of the two. Also, it can be noted in figure 8 that the relative 
advantage of one data type over the other is not constant, that is, the curves 
diverge. It is concluded that the data type producing the best set of elements 
is very dependent on the semimajor axis. 

A range of satellite orbital inclinations from 2^ to 40^ with respect to 
the earth -moon plane was considered. The elements a, n, 03, e, and to were 
held constant at the nominal values while i was varied over the given range. 

It was assumed that the satellite had been tracked over a period of five con- 
secutive orbits by making 26 range and 26 range-rate observations per orbit 
equally spaced in time. The results are presented in figure 9* This figure 
shows typical curves for the variation of the one -sigma error in determining 
the elements with the sine of the inclination. Over the range of inclinations 
in the investigation, the accuracy of determining the orientation angles i, 

Qy and O) increased significantly as the inclination increased, whereas the 
accuracy of determining the in-plane elements a, e, and to remained approx- 
imately constant. It should be noted that the curves for and have a 

slope of -2 and the curve for has a slope of - 1 . Similar plots not pre- 

sented herein showed that these slopes are independent of the nodal position 
and hence it was concluded that Gq and are inversely proportional to 

sin% and that aj_ is inversely proportional to sin i. This result was 
unexpected inasmuch as it was not apparent from the form of the partial deriv- 
atives. It can be seen from figure 9 that over the range of inclinations from 
2 ^ to 4o^ either range or range -rate measurements can be used to determine the 
elements with approximately equal accuracy when ^p/o'p = 1500 seconds. It was 
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concluded from this figure that the relative advantage of one of the data types 
over the other does not change over the range of inclinations in the present 
investigation . 

The curves shown in figure 10 indicate that the position of the line of 
nodes during the tracking period is a very significant parameter in the deter- 
mination of the orhital elements. The results presented in this figure are 
those obtained by assuming that the satellite had been tracked over a period 
of five consecutive orbits during which time the nodal line had rotated through 
an angle of 7^. The values of Q, shown in figure 10 are the values at the 
beginning of the tracking period. Again 26 range and 26 range-rate observations 
were assumed to have been made during each orbit ^ with the elements other than 
being held at the nominal values. For the elliptic orbit in the present 
investigation^ it was found that the elements 0.^ cd, and e are best deter- 
mined when the orbit is viewed on edge (n = 0^ and l 80 ^) whereas the elements 
a, i^ and are best determined when the orbit is viewed broadside 

(^7 = 90 ^ 270 ^). The very large variation in the curve for is due to 

the fact that^ as stated earlier^ the partial derivatives of range and range 
rate with respect to i are approximately proportional to sin 0. As Q 
approaches 90 ^ stnd 270 ^^ the amplitudes of these derivatives approach the maxi- 
mum values^ and therefore the inclination is more accurately determined at 
these nodal positions. The curves in figure 10 show a periodic property with a 
period equal to 1 lunar month. They are also symmetric about 90^ a,nd l80^ 
because occultation was not considered. If occultation had been considered^ 
this symmetry would have been partially destroyed. In the discussion of fig- 
ures 2 and 5 it was stated that bp/hQ is approximately equal to bp/bod and 
that bp/bO, is approximately equal to 8p/6cD for near-zero inclinations. 
Therefore for small inclinations, high correlations between 0 and cd and a 
similar accuracy of estimation would be expected. It can be seen in figure 10 
that even with a mediiom inclination, the curves for Cq and are very much 

alike with approximately the same variations. The correlation coefficient 
after five orbits is 0 . 9985 * It can also be seen in figure 10 that the relative 
advantage of one data type over the other is not constant over the range of 
nodal positions in the investigation but remains within a factor of approxi- 
mately 2 . 

The effect of eccentricity upon the accuracy of estimating the elements 
is illustrated in figure 11. It was assumed that the elements i, H, o), and 
to were held constant at the nominal values while the semimajor axis was fixed 
at 5000 kilometers and e was varied between 0.01 and 0.6. The results shown 
are those obtained by ass-uming that the satellite was tracked over a period of 
five consecutive orbits with 26 range and 26 range-rate observations made per 
orbit. Except for the accuracy of estimating e which had only a slight var- 
iation, the accuracy of estimating the elements was found to increase as the 
eccentricity increased. However, it can be seen that for eccentricities above 
0.1 the acciiracy of determining the orientation angles i, H, and co is not 
appreciably improved when e increases whereas the scale of the orbit, which 
is inferred by the semimajor axis a, improves significantly. The results 
indicate that the accuracy of estimating the time of periapsis passage to is 
very dependent upon the eccentricity. Since the position in orbit is dependent 
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upon and^ furthermore^ since the knowledge of to greatly improves with 

an increase in eccentricity^ it can he concluded that the position in orhit is 
better determined at the higher eccentricities. It should also he noted that 
the data presented in figure 11 seem to indicate that range measurements are 
better than range-rate measurements for all values of e; however^ it must be 
remembered that in this case the semimajor axis is fixed at 5000 kilometers and 
that, as noted previously, range data have a definite advantage over range-rate 
data for large values of a. Hence, it should not be concluded that range meas- 
urements are always better than range-rate measurements for all values of e. 
However, it can be concluded that for the elements a and to the relative 
advantage of one data type over the other is somewhat dependent upon e. As 
shown in the following paragraph, the accuracy of determining the position of 
a satellite in orbit is primarily limited by the accuracy of determining the 
orientation angles, and it can be seen from figure 11 that the relative advan- 
tage of one data type over the other for determining these angles is independ- 
ent of e. Therefore it can also be concluded that in the determination of the 
position of a satellite, the relative advantage of one data type over the other 
is independent of e. 

One of the more important parameters to be estimated in any orbit- 
determination problem is the position of the satellite. It is therefore of 
interest to know how accurately this position can be determined. A qualitative 
expression for the accuracy of determining the position in terms of the one- 
sigma errors in the estimates of the elements can be obtained as follows. The 
vector r is written as 


r = r(a,i ,f2,U),e,to^t) 


( 10 ) 


and therefore at some fixed time 


A? = ^daH-^Ai+^An + ^^ + ^Ae + ^Ato (ll) 

Sa Si Sq Se St q 

Equation (ll) is then put into rectangular (x,y,z) component form where the 
partial derivatives of x, y, and z with respect to the orbital elements are 
maximized with respect to the position in orbit and the angular variables i, 

Qy and cjd. By using the triangle Inequality, the length of the vector can 

be written as 


1^1 < |Ax| + |Ay I + 1 Az| (12) 

The maximized partial derivatives can then be substituted into the relation 
given in inequality (l2) and if the largest coefficient of each incremental 
change in the elements from the three-component inequality is selected, the 
following equation can be written: 
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I 


1^1 ^ 5 


(15) 


(l + e)Aa + a(l + e)(Ai + ZSQ + 


+ 



(3 Ae + n Ato) 


By assuming that the small incremental changes in the elements are equal to the 

stan(iard deviations of the elements and that |z^j is equal to the standard 
deviation of the position F, inequality ( 15 ) can be written as 


r 


CTp 


<3 


(1 + e)cTa t a(l + e)(cTi + + o^) 



( 1 ^) 


It can be shown by substituting the one-sigma errors in the estimates of the 
elements from figure k into inequality (l4) that the major contribution to the 
uncertainty in the position is the uncertainty in the orientation angles. 

Hence ^ it is very important to make measurements which allow an accurate deter- 
mination of these angles in order that the position of the satellite might be 
detennined with the most accuracy. 

Another point of interest, which is also noted in reference 3^ is that the 
reflection of any given orbit through the earth-moon plane would give the same 
time history for range and range-rate measurements as that of the original 
orbit. This fact implies that an addition of l80^ to both the node and the 
argument of periapsis would not result in any change in the accuracy of deter- 
mining the elements and that without a priori information it would not be known 
which of the two orbits was being tracked. 


CONCLUDING REMAEKS 


A parametric stu(3y of the effects of tracking mode and orbital parameters 
on the accuracy of determining the state of a lunar satellite has been made by 
using range and range-rate measurements. A comparison of the relative advantage 
of using either range or range-rate measurements indicates that the advantage 
of one data type over the other is very dependent on the semimajor axis and 
somewhat dependent on the nodal position and eccentricity but is not dependent 
on the inclination. For Ixmar orbits with meciium eccentricity and with semi- 
major axes of approximately 25 OO kilometers, the orbital elements could be 
deteimined with equal accuracy by using either range or range-rate data when 
the ratio of the standard deviation of range measurements to the standard 
deviation of range -rate measurements was approximately I 5 OO seconds. For lunar 
orbits with semimajor axis smaller than approximately 2500 kilometers, range- 
rate measTxrements gave a more accurate determination of the elements, whereas 
for orbits with semimajor axes greater than approximately 2500 kilometers, range 
measurements proved to be the better data type of the two. It was concluded 
that range and range -rate measurements are similar data types in that they pro- 
duce similar correlation matrices for the state variables and simultaneous 


Ik 



use of both data types does not significantly reduce correlations between 
the elements. 


Langley Research Center, 

National Aeronautics and Space Administration, 
Langley Station, Hampton, Va., October 11, 
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APPENDIX A 


DETAILS OF WEIGHTED LEAST SQUARES SIMULATION 


A special case of a weighted least squares process was used in the error 
analysis presented in this report. This appendix is included to illustrate 
the basic equations which were used in the IBM 709^ electronic data processing 
system. In particular, the method of obtaining the covariance matrix of the 
orbital elements is shown. This method is essentially the same as that 
described in reference 4. 

The solution of the equations of motion of a point mass about a central 
body contains six constants of the motion which may be taken as a, i, cd, 

e, and to* The observable quantities, range and range rate, written as func- 
tions of these six constants at any given time t are 


P 

P 


f(a,i,a,CD,e,to^t)| 

F(a,i,a,cu,e,to,t)| 


(Al) 


Theoretically only six properly chosen measurements would be required to deter- 
mine the six elements of equations (Al) provided there were no errors in the 
measurements. Since a measurement error is associated with any measure- 
ment or more than six measurements can be used to obtain a "best" set 

of elements. The notion of "best" is to be defined subsequently. The sub- 
script i denotes the ith measurement. 

Equations (Al) represent the functional equations for range and range rate 
but similar representation can be made for any measurable quantity. The fol- 
lowing equation can be written for the ith measurement of the general quan- 
tity y and error in the ith measurement e ± : 



In eqiiation (A2), is analogous to the ith measurement of range or range 

rate (defined in eq. (Al)) whereas a. is analogous to the elements in equa- 

tion (A1) . In general, f is a nonlinear function in the a. 's, and in order 

to make the problem amenable to solution, the basic equations are linearized 
about a nominal set a9. If the true values of aj are assumed to be close to 

a nominal set a9, then y^(t^) = y? = f(cc9,tj_) where y? is the calculated 

J Id -L 

value of the ith measurement obtained by using the nominal set aj . The fol- 
lowing equation can then be written; 
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A/^i = - y? = + ei (AJ) 

If eq^aation (AJ) is expanded in a Taylor series about the nominal set aj and 

if the terms of order higher than the first are dropped^ the following equation 
is obtained: 






(a4) 


The notation 


J 


means the partial derivative of 



with respect 


lo 

to the jth element evaluated at the prescribed set a?, j = 1^2^ . . .6^ at 

o 

the time ti . By denoting = Ax j , equation (a 4) can be written as 




Sf 


(ti) 


hcLA 

0=1 ^ 


AXj + €i 


(A5) 


Equation (A5) can be written in matrix form as 


Ay = A Ax + e 


(A6) 


where Ay is an N x 1 column vector of the observed minus the computed val- 
ues of Yy A is an N x 6 matrix of known partial derivatives^ Ax is a 
6x1 column vector of deviations of the elements from the nominal set^ and € 

is an N X 1 coliimn vector of observation errors. The problem is to find the 
A A 

best estimate Ax of Act. when Ay and A are given. If Acc is determined 

and the nominal values aP are used^ the best estimate of a is a = + Aa.. 

In order to specify what is meant by best estimate^ some quantities must 
first be defined. Denote a as the best estimate for cLy y = Aa as the best 
estimate of the true value of the observable, and y - y = y - Aa = e as the 
best estimate of the observation error. Note that the deltas have been dropped 
to simplify notation. 

The best estimate of a is now defined as that a which minimizes the 
sum of the squares of the weighted components of the residual vector y - Aa. 
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In order to account for the difference in confidence between various observa- 
tions and the possible relations between them^ a so-called weighting matrix W 
which is assumed to be a symmetric^ positive -definite N x N matrix is intro- 
duced. Ihe wei^ted least squares function which is to be minimized can now 
be written as 


Q(o') = (y - Aa)%(y - Aa) = 


(AT) 


In order to minimize the function in equation (A 7 ), the variational prin- 
ciple given in reference k is used. This principle states that in order for 
Q to be an extremum^ the first variation in Q must vanish and in order for 
this extremum to be a minimum^ the second variation must be positive-definite. 
If this principle is applied^ with only the variation in oo being considered^ 
the following equation can be written: 


6Q = -6a^ A%(y - Aa) - (y - Aa)%A5a 
= -26a'^ A%(y - Aa) 


(A8) 


The value of this equation must be zero for an extremum^ and the fact that 6a 
is arbitrary implies that the best estimate a must satisfy the equation 


A^W(y - Aoc) = 0 


(A9) 


or 


A%Aa = A^y 


(AlO) 


Premultiplying both sides of equation (AlO) by gives the best 

estimate 




a = 



(All) 


provided a"^WA Is nonsingular. In order to show that this is the a which 
minimizes it is sufficient to show that the second variation is positive- 
definite^ where the second variation is 


62q = 26a^ A^WASa 


(A12) 


For arbitrary Sa, equation (A12) is greater than zero if W is positive- 
definite. One of the basic assumptions was that W is positive-definite^ and 
therefore a in equation (All) is the best estimate of a. 
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The remainder of this appendix is devoted to the development of the covar- 
iance matrix. Equation (All) has been used to determine the best estimate a 
and it is now desirable to determine the statistics of a. For example^ how 
well was a estimated? The definitions of the variance and the first mixed 
moment can be used to write that the covariance matrix for a is ECoa?*) 
where E denotes the expected value of the variable in the parentheses. By- 
use of equation (A11)^ the following equation can be written: 


E 


= E (a'I\7a) \%yy%'^A(A%^A)~^ 


(A13) 


Since all terms except the random variable y are constant in equation (A 15 )^ 
the operator E operates only on y and yT. Hence ^ 


E(oa^) = (a’^Wa) ^A^^^WE^yy'^jw^A^A*^ (A1^) 

Reference k shows that the best choice of W is the inverse of the covariance 
matrix for the measurements where this covariance matrix is E(yy^). Assume 

that each measurement is of equal weight . e . , . Oy^ and 

that the observations or measurements are completely uncorrelated. If these 
assumptions are made and if it is noted that E(yyT) is the covariance matrix 
for the measurements^ then 


W = 



(A15) 


and 


WE(yy'^) = I 


(Al6) 


where I is the identity matrix. Hence^ equation (A14) becomes 

= E(aa.T) = (a'IVa)"4a%^a)(a%'^a)"^ 

= (aTwa)-I = a/(ATA)-l 


(AI 7 ) 
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Note that in the terms of the orbital elements 


AS. = 


cJa^ 



PaCD^a^o) 

PaePa^e 

Pato^a^to 

Pla°i®a 




PlePi^e 

Pito^'i'^to 

^na'^n'^a 










PcuePajPe 

^a^Q^oPto 

•^ea'^e^a 

Pei°ePl 



Pe^ 

PetiQ^s^'to 

Ptoa°to°a 




PtoePto°e 

04 

0 



(Al8) 


where Pai^a^i-^ ^an^a^n<? forth; are the first mixed moments and Pal.? 

pg^^; and so forth; are the correlation coefficients. This matrix was used to 

obtain the results in the present report. It can be seen from equation (AI7) 
that no actual measurements are needed to determine the variances of the 
elements; in particular; only the variance of the measurements is needed. The 
covariance matrix for the simultaneous use of two data types such as range and 
range rate can be written in several wayS; one of which is 


a 




-1 


The correlation matrix is a matrix having ones as the diagonal terms and 
the correlation coefficients of equations (A 18 ) as the off-diagonal terms. 

It should be remembered that equation (AI7) is a special case of the gen- 
eralized weighted least squares as given in equation (AI3) . From the special 
form given in equation (AI7); it can be seen that for a single data type the 
variances of the elements are proportional to the variance of the measurements. 
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EQUATIONS FOR PARTIAL DERIVATIVES OF RANGE AND RANGE RATE 
WITH RESPECT TO ORBITAL ELEMENTS 


The elements of the A matrix which are the partial derivatives of range 
and range rate with respect to the orbital elements are as follows: 


I - i((r * 


n(t - to) sin E 
2r 




Rl2 5 n(t 


to)(l - e2) 


1/2' 


— = i Rrni sin Q 

bi P ^ 


(Bl) 


(B2) 


^ = . 1 Rrmi 


(B3) 


. i Er!2 

5o) P 


(Blf) 




sin^E - cos ^ 


RrZ2 




sin E 


(1 - e2) 


1/2 


+ (1 - e2) 


^/^sin ^ 


(B5) 


bt, 


~ ® ^^2 ^(1 - e 2 )l /2 


(b6) 
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where 


p2 = r 2 + r2 + 2RrZi 


P 


I ^ 

Sn Sto 


^ — n(t - to) sin E 

aa a 2 ^ 

^(E - M) = - 2 n(t - to)- cos E 

aa 2 I* 

^ = - I ^ ^ 2 n(t - to)(l - e 2 )l /2 

^ = - I ® 2 n(t - to)(l - e2)l/2 

c3sl 

^ = I ^in(t - to)(l - e 2 )l /2 
oa 

In the equations in this appendix a^ i^ co, e^ and to are the conven- 

tional Keplerian elements. The mean anomaly M is given hy Kepler’s equation 

M = n(t - to) = E - e sin E 
where n is the mean angular rate of the satellite. 

The direction cosines are as follows: 

= cos 0 cos Q - sin 0 sin Q cos i 
l2 = -sin 0 cos Q - cos 0 sin H cos i 
mj^ = cos 0 sin Q + sin 0 cos 0 cos i 
m 2 = -sin 0 sin ^2 + cos 0 cos Q cos i 
nj_ = sin 0 sin i 
n 2 = cos 0 sin i 
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Partiol derivatives of range with respect to the elements 
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Figure 5.- Variation of the standard deviation of the nominal-orbit elements with the number of observations. The observations were assumed to 

have been made over five orbits. 
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Figure 9.- Variation of the standard deviation of the nominal-orbit elements with the sine of the inclination. 
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Figure 11.- Variation of the standard deviation of the nominal-orbit elements with eccentricity. (Note that a was fixed at 5000 km.) 
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